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We consider a driven single mode Dicke-Hamiltonian coupled to a dissipative zero-temperature 
bath. We derive the cumulant generating function for emitted photons of this quantum-critical 
system by using a P-representation of the master equation in the thermodynamic limit. This 
cumulant generating function is shown to consist of two parts: a macroscopic component, which is 
Poissonian in nature with characteristic rate proportional to the order parameter of the system; and 
a part describing fluctuations which is non-trivial in form and divergent around the quantum phase 
transition. 
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PACS numbers: 42.50Ar, 42.50Lc, 05.30 Rt 



I. INTRODUCTION 



The Dicke model describes the interaction of N two- 
level systems with bosonic field modes. In the thermo- 
dynamic limit, and for a single bosonic mode, it un- 
dergoes a quantum phase transition (QPT) from a nor- 
mal to a superradiant phase when the atom-field cou- 
pling strength exceeds a critical valu^^HH, xhe model 
has been experimentally realised by an B ose- Einstein- 
condensate trapped in an optical cavitj^^. Dissi- 
pation due to photon emission has furthermore been 
modeled, e.g. by Heisenberg-Langevin equation^^ltiSl^ 
the Keldysh approach?^ or by use of Hartree-Fock- 
Bogoliubov theorjEl 



In this paper, we study the full photon counting statis- 
tics of the driven dissipative single-mode Dicke model 
by including a counting field x in the master equation, 
which we analytically solve in the P-representation in 
the thermodynamic limit N oo. We obtain analytic 
expressions for time-dependent cumulants of the photon 
counting statistics, as well as the asymptotic cumulant 
generating function (CGF). The CGF is shown to consist 
of two parts: a macroscopic component, which is Poisso- 
nian in nature and has a characteristic rate proportional 
to the mean occupation of the cavity mode (order param- 
eter of the system); on top of this comes a contribution 
describing fluctuations about the mean behaviour, which 
has a non-trivial form and is divergent around the quan- 
tum phase transition. In addition, we identify the three 
phases (normal, superradiant, intermittent) and the cor- 
responding critical coupling parameters Ai,A2,A3 that 
characterize the dissipative phase transition. 



The structure of this paper is as follows: In Sec. [Tl] we 
introduce the model and obtain the three critical coupling 



parameters, in Sec. Ill we solve the master equation with 



a counting field Xj in Sec. IV we give an analytical ex- 
pression for the CGF and in Sec.|V]we discuss our results 
and connect them to other works. 



II. MODEL 



The Dicke Hamiltonian {h = 1) 



A 



H — ujqJz + ujolai ^ — (al + ai)( J+ + J_) (1) 

describes our isolated system. Here, , J± are the collec- 
tive atomic angular momentum operators that describe 
an ensemble of N two level atoms with a level splitting 
ujq; j is the length of pseudo-spin with value j = N/2 for 
bosonic realisatiorP^, oi and a\ are the ladder operators 
for the optical mode with an energy w, and A is the cou- 
pling strength between the optical mode and the atoms. 
As in Refs.^ we interpret Eq. ([T]) as an effective model in 
a frame rotating at the frequency of an external driving 
laser, ^l. 

We describe dissipation in this system with the master 
equation in Lindblad-fornP 



dt 



p{t)^~i[H,p{t)] 



i\aip{t) — 2aip{t)a\ + p{t)a\ai 



(2) 



where p is the density matrix of the atom-cavity system 
and r is the rate of photon loss from the cavity's! 
To make analytical progress, we use the Holstein- 
PrimakofF-transformatiortl^, which allows one to describe 
iV-two-level-atoms by a single bosonic mode. This is done 
by the transformation 



J+ = 4\/2j - 4^2; 
Jz = (4^2 - j), 



(3) 



where 02 and a\ ladder operators of the intro- 

duced atomic mode. In this representation, the Dicke- 
Hamiltonian reads 



H = ujo{a\a2 — j) + wajai 



+\{al + fli) a 




4«2 



4«2 



«2 .(4) 
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A. Thermodynamic Limit N ^ oo 

The isolated Dicke-Hamiltonian H in the thermody- 
namic hmit iV — >■ oo is known to have two phases: normal 
and superradianlP. 

In the normal phase and for N oo only H becomes 
effectiveljB', 

iJjv = '^oO-2'^2 + LiJa\ai + X{a\ + 02) (a| + ai), (5) 

which has the form of two simply-coupled harmonic os- 
cillators. The master equation for this phase may be ob- 
tained by simply replacing H in Eq. ^ with this form. 
In contrast, the superradiant phase is characterised by 
the macroscopic occupation of the ai and a2 modes. In 
order to describe this phase, we insert a mean-field-ansatz 

fli = c+\/a, a[ = c'+^/a*, 02 = d— 02 = d^ — \/p*, 

(6) 

into both H and the dissipator term of the master equa- 
tion, Eq. ([2]). Let us start with H, where we insert Eq. 
^ into Eq. Q. We then expand the square roots and 
neglect terms with powers of N in the denominator'*. As- 
suming G IR we obtain in the thermodynamic limit 



Hs - Hg 



(1) 



H 



(2) 



(7) 



with 




rZ 



(2) 



2j 

LOC' C + 

+Mid^ + df, 



(8) 
(9) 



with 



k = 

A = 

M = 




2fc2 



(10) 



Making the displacements in the dissipator, the com- 
plete master equation then reads 



dt 



pit) = -^ 



ac^ + i — ^/a*c^ p{t) 

{c^cp[t) + p{t)c^c ~ 2cp{t)c^) 

(11) 



.r 



Now, we first determine the macroscopic occupation 
\a\ , 1/3 1 in the superradiant phase using the master equa- 
tion. In Eq. (11) we have already separated the terms 
that are linear in the new operators c, . These terms 
are proportional to the square root of the particle number 
N, so they diverge in the thermodynamic limit. We de- 
termine the macroscopic occupation, such that all parts 
~ \/N of the master equation vanish. That means the 
first commutator has to be zero which is fulfilled if 



H 



(1) 



i — y/ac' 



' a c ■ 



0. 



(12) 



Next, we insert the H^"^ from Eq. (|8]), factor out the 
operators c, , d, d^ and get an expression of the form 
c(..) -I- c'i..) + d{..) + d\..) = 0. To achieve the identity 
we set all braces to zero and arrive at four equations. 
The last two of them are equal, because we have assumed 
/3 e R. These equations can be solved (see Appendix [b|) 
to yield 



2A- W2j 1- 



• A2 

, A 



13 = J I 




itlPl 



(13) 



(14) 



which differs from the isolated case by the inclusion of 
the loss-rate F. For A > A2, a and /3 are non-zero. The 
coupling A2 is our first of three critical points. It is the 
coupling above which the macroscopic excitation exists. 
There is also one trivial solution, a = /3 = 0, that recov- 
ers the normal phase. 



B. Analysis of the non-dissipative part of the 
master equation 

Some insight into the behaviour of the system can be 

obtained by just considering the non-dissipative part, i.e. 

(2) 

the Hamiltonians Hjq and Hg with displacements cho- 
sen as in the foregoing. We study here the eigenvalues of 
these two Hamiltonians and in the next section we will 
rewrite the master equation, Eq. ([2]), in corresponding 
diagonal bases. 

Both Hamiltonians can be diagonalised by a transfor- 
mation of the following forrcP 



/ A 


B 


G 








B 


A 


D 


a 




d\ 


M 


B2 


G2 


D2 




d2 


\B2 


A2 


D2 


G2 j 




\4J 



where 



(c, c\d, d^)^, 



normal phase, 
superradiant phase 



(15) 



(16) 
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represents the old basis and di, the new basis, in which 
the system is diagonal. In both cases, the diagonalis- 
ing matrix has the same structure but different internal 
parameters, which are listed in Appendix Eq. ( Al )-( A6 ). 

(2) 

After diagonalisation, Hs or Hg assume the form 



Hd — e-d\di + £j^d\d2 with eigenenergies 



p(i) 

^(2) 



in the normal phase, 

in the superradiang phase 



the forms of which are given in Eqs. (A2),(A4| 



(17) 



In the normal phase, one of the eigenenergies e"^^ of 
Hjsi has a zero at A = Ai with 



2 _ UJUJq 



(18) 



The vanishing eigenenergy marks the end of the normal 
phas^. 

Coming from the superradiant phase, the same 
eigenenergy vanishes, but in this case at a critical value 
of A = A3: 



A^ 



+ 4^2) Wo 

16 • 2^2 



A? 



+ 0{T^). (19) 



We recognize that A3 and the eigenenergies in the su- 
perradiant phase also depend on the bath-coupling- 
parameter r. Furthermore Ai < A3, so we obtain a gap, 
where the eigenenergies become complex and our analy- 
sis breaks down, see the discussion below. In Fig. [Tj the 
excitation energy of the closed system (case F = 0) and 
open system (F 7^ 0) are plotted. 



III. MASTER EQUATION IN 
P-REPRESENTATION WITH COUNTING FIELD 

In the following, our aim is to count how many pho- 
tons are lost to the environment, and to determine the 
properties of the photon distribution function. To this 
end, we introduce a counting fielcP^ x i^ito the master 
equation ([2| 



dt 



P{x,t) = -i[HN/s,p{x,t)] 



a{aip{x, t) - 2e''^aip(x, t)a\ + p{x, t)a\ai 



(20) 



To make our calculations analytically tractable, we make 
a rotating wave approximation (RWA) . We transform the 
master equation, Eq. (20), to a diagonal basis di of 



(2) 

in the normal phase or of Hg in the superradiant phase 
using the relations 



ai = Adi + Bd\ + Gd2 + Dd] + ^/a , 
a\ = Bdi + Ad\ + Dd2 + Gdl + Va* , 



(21) 
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FIG. 1. (Color online) Eigenenergies of the effective non- 
equilibrium Dicke-Hamiltonian from the master equation 
plotted for both phases as a function of A/A2 (solid lines). 
Also plotted are the energies for the equilibrium (F = 0) case 
(dotted lines). For A < Ai we have a normal phase, and for 
A > A3 a superradiant phase. At A = Ai and A = A3 the e_ 
energy vanishes. For F = 0, Ai and A3 coincide, but out of 
equilibrium, a gap appears. Parameter: = 2F,a; — F (case 
F = 0: tJo = 2uj). 



where the coefficients are defined in Appendix [Xj Next, 
we transform into the interaction picture and neglect all 
fast rotating terms such as did^j for i j. We obtain the 
following RWA-master equation for p = p(x, t) 

F 



A^idldip - 2e'^dipdl + pd\di) 

B^{did\p - 2e'^d\pdi + pdid{) 
D^{d2dlp - 2e^^d\pd2 + pd2dl) 

G^{d\d2P - 2e'^d2pd\ -V pd\d2) 



F|a|(l 



Op- 



(22) 



This master equation describes two uncoupled harmonic 
oscillators in interaction with two independent thermal 
baths at different temperature^^. We see that each ex- 
change of quanta with these baths, either into or out, is 
associated with the emission a of 'physical' photon (i.e. 
one with annihilation operator a\). The final term in 
Eq. ( 22 ) is only nonzero in the superradiant phase and 



represents a Poissonian process with rate proportional to 
a macroscopic excitation \o\. 



Now we transform the master equation ( 22 ) into the 
P-representation using the following ansatz 

P = / f'(7i,7i,72,72:X,^) I71) (71I «) I72) (72I ^^710^^72 , 



(23) 

where [7^) is the eigenstate of the annihilating operator 
di with eigenvalue 7i. Using diYi^ = ,d\\^i) = 



4 



{-^ + 7*) |7i) and integration by parts, we get a par- 
tial differential equation for the quasi distribution P = 

P = 6P-T\a\{l-e'^)P, (24) 



However, we repeated our calculations (x = 0) in 
the singular coupling limit of the master equation where 
terms like e"*"*^** approach unity and found that, for the 
parameters studied, no significant qualitative differences 
with the RWA method arose. 



with 



o = u 



d'fi 



7i + ^7rj+v 



-72 



:l2 



r(i 



{{A^ + 52)717* + (G^ + 52)727* + 



52)1; 



2e^^)) ; 



Equation (24) can be solved as 

P = Fo(t).e-i^l"l(i-^'")*^ 

if 

Po(i) = OPoit) 



(25) 



(26) 



(27) 



is fulfilled. Since O is a bi-linear operator, the solution 
of Eq. (27) can be written in the following form 



Po = cxp(-a(x,t) + bi(x,t) ■ 71 + 62(x,t) • 72 
+ ci{x,'t) -71* +C2{x,t) -1*2 
- di{x,t) ■ 7i7i* - d2{x-,t) ■ 7272) • 



(28) 



Substituting (28) into Pq ~ OPa wc obtain seven cou- 



pled first order differential equations for the functions 
a,bi,Ci,di, see Appendix Eq. (B8). Taking a displaced 



2-dimensional gaussian distribution with a standard de- 
viation £ as an initial condition. 



-fo(7i,7i,72,72>X,* = 0) 



1 



(27r£)2 



(29) 

we can solve this system analytically. The steady state 
solution for t — > 00 can also be calculated and is listed in 



Appendix Eq. ( B9 ) 



With knowledge of the P-representation of the density 
matrix with the counting field x, we can now calculate the 
mean occupation of the original modes and the cumulant 
generating function for the photon counting statistics. 

As a final remark, the RWA would normally not be 
valid exactly near the points Ai,A2, because here one of 
the eigenenergies e± vanishes. Terms like Ji^^'di^'' can no 
longer be neglected on the basis that they are quickly 
rotating. 



IV. RESULTS 



A. System properties 



We first calculate the occupation number of the pho- 
tonic and atomic modes. In our P-representation, an 
operator average is given by 



TvAp = Tr / AP I71) (71I (g, I72) (72I djfdji 



(30) 

which we use to calculate (a\ai^ and ^0202^ We have 
to use the relation between the old and the diagonal 
basis, because the P-representation has been evaluated 
in the diagonal basis. We also have to do the calcu- 
lation in the interaction picture and use the RWA, as 
we have done in the calculation of the P-function. In 
both phases the time-dependent average of the optical 
and atomic modes has the same structure of the form 



00) from which we can read off the characteristic relax- 
ation times 



n = (rcos2 7)-i and T2 = (rsin2 7)-i , (31) 

where 7 is the rotation angle of the decoupling, Eqs. 
(|A3|),(A6). As a result, the two modes develop at two 



different time scales. In the following, we look at the 
long-time solution. 

In the normal phase, using the long-time solution of 
Eq. ( 24 ) we get (see Appendix [C]) 



"^"0^'^°")=2c...o-8A2' 
(4a2)(t^oo) = i(-4+^ + 



2cjo 



(32) 



from which we recognize the expected divergence at A = 
Ai, Eq. (flSl). 



In the superradiant phase, we have a finite macroscopic 
occupation, which is much greater that the fiuctuations 
around it. We obtain 



a\ai^ = (c^c) + \a\ ; 



(33) 



5 



where in the stationary case t ^ oo 



io"^ sin^ 



(7) 
2 



=(2) 



(34) 



w ilo wriosin^(7) 



COS' 



'ii) 



+ 



(2) 



=(2) 



ujHq cos^(7) 



(35) 



Fig. [2] shows the occupation number of the field and 
the atomic modes as a function of A in the steady state 
and gives results that are approximately equivalent to 
those of Nagy and co-workers^. We see that the fluctu- 
ations diverge around the phase transition, but there is 
again the undefined area for A € [Ai, A3]. This area has 
been also seen in the energy plot. 

In the normal phase the occupation of both modes de- 
crease with smaller coupling parameter A. Decreasing A 
weakens the coupling between the optical and the field 
mode. For A — ^ the atomic mode decouples completely 
and only the field mode becomes damped. 

In the superradiant phase, the occupation of both 
modes first decreases with increasing A but then the oc- 
cupation of the field mode increases. The increasing A 
reduces the coupling A, Eq. (10 1, which explains the de- 
creasing. On the other hand, the A/-term in the Hamil- 
tonian then becomes dominant. For A A3 the d-mode 
becomes approximately decoupled but squeezed because 
of the presence of the M-term, see Eq. Q, which ex- 
plains the increase of occupation. That means, only the 
c-mode is then damped by the bath. 

The macroscopic occupation \a\ , |/3| exists only in the 
superradiant phase, where it is much greater than the 
fluctuation, and vanishes for A < A2. 



B. Counting statistics 

We now derive the cumulant generating function 
(CGF) F{x,t) of the photon counting statistic^, de- 
fined as 



F(x,t) =logTrp(x,t). 
The fc-th cumulant can be obtained by 



(36) 
(37) 
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FIG. 2. (Color online) Fluctuation of the field and atomic 
mode occupation number (solid line), Eq. (32 1,(33 1 and the 
correspondent macroscopic occupation (dashed line). User 
parameter: uj = 2r,a;o = O.SF. 



In the long-time limit, the CGF can be calculated from 
our P-function including a counting field x, Eq. ( 28 ) as 
follows 



F{x,t^ 00) = log J P{oi, 0*1,02, 0*2, Xit (X))(foi(fo2 
= -T\a\[l-e''^)t 



log 



4^^ 



^di{x,t)d2{x,t) 
~T\a\{l~e'^)t~a{x,t^oo), (38) 



where a{t 00, x) is defined in Appendix Eq. (B9) 
and terms di <^ a have being neglected. Inserting the 
quantity a(x, t 00) from Eq. ( B9 1, we obtain the CGF 
for both phases, 

lim -F(x, t) = T \a\ (e'^ - I) + U - + & - 

i— T-oo t 2 



A-^ + m- 2^252 (2e2*x _ 1) 

L>4 + (74 _ 2D'^& (2e2*x - I) ) F . 

(39) 



The CGF thus consists of two parts. The first is the 
CGF of a Poissonian process with rate F|a|. In the 
superradiant phase, \a\ is macroscopic (proportional to 
the number of atoms) and this contribution dominates. 
The k-th cumulants from this macroscopic contribution 
are simply (n^^ ^ = T \a\. On top of this contribution, 
which is absent in the normal phase (where a = 0), 
comes a further contribution to the CGF which arises 
from fluctuations about the mean-field displacements. 
The cumulants due to this contribution are given by 

With this result we are now able to calculate all cu- 
mulants in the long-time limit. For example, the first 



FIG. 3. (Color online) Asymptotic cumulants divided by time 
in both phases for driven dissipative Dicke model, {n^) ^ ~ -'V 
is macroscopic, 5 (?^'°)^ is the corresponding fluctuation, Eq. 
( 40 1 . For A 2> A2 the cumulants starting by the third are 
growing. We can see also a tiny gap around A2 corresponding 
to the undefined area Ai < A3 of our theory. Parameters: 
uj = 2V,u)Q = o.5r. 



FIG. 5. (Color online) Temporal development of the first 
five cumulants and of the photon mode occupation average 
(^a\a-i^. After some characteristic time scale, (^a\a\^ ap- 
proximately reaches the steady state value, and on the same 
time scale the cumulants change their slope. Parameters: 

ujo = 2r,u = ir,A = 0.3/r. 



Fano -Factor 




FIG. 4. (Color online) Same as Fig. [sjbut with Fano factors 
not cumulants. 



cumulant reads 



(">c = |a|rt + (5(n)p , 



where 



ic 



= Tt 




normal phase 
superradiant phase. 



(40) 



(41) 



We see that 5 {n)(j is proportional to the mean occupa- 
tion of the optical mode that directly couples to the envi- 
ronment. This fits with the experiment use of the mean 
macroscopic flux ^ |a| as a measure of the o rder parame- 
ter of the superradiant phase transitiorPsSI, As must be, 
the cumulant generatin g fu nction is linear in time. We 
mention that the source^ of the photon current is the 
driving frequency fl. 

In Fig. |3] we have plotted the first five cumulants di- 
vided by time as a function of the coupling parameter A. 



We see all cumulants diverging near the phase transition 
A « A2, where one eigenvalue of the Hamiltonian be- 
comes zero. Going away from the critical region, our sys- 
tem emits less to the environment and therefore 5 (^n'^) /t 
decreases. For A A3, J (n)^ , S {n'^)^^ become constant, 
however, higher cumulants such as (5(n^)^ grow. How- 
ever, in the superradiant phase all fluctuations can be 
neglected because of the existence of a macroscopic Pois- 
son distribution {n'')^ (Note the scaling 1/N in Fig. |3|, 
which is much greater than the fluctuation. We see again 
the undefined area around A2 . Fig. [4] shows the fluctua- 

tion part of the Fano-Factors Fi = ■ 

If we use the complete solution of Eq. ( p4| ), we can 
calculate the cumulants as a function of time t. The 
evolution of the first cumulants is shown in Fig. [5] . We 
see that there are at least two time scales. This behaviour 
is due to the two timescales ti and T2, Eq. ( [sTj ), of the 

both modes. As we can see in the same picture, (a\ai 
develops in an analogous fashion. 



V. DISCUSSION 

We have derived the cumulant generating function for a 
driven single-mode Dicke Hamiltonian coupled to a bath. 
The combination of the thermodynamic limit and the use 
of the P representation allows us to obtain exact expres- 
sions not only for the asymptotic CGF and cumulants, 
but also those at finite times. 

The CGF and photon-counting cumulants consist of 
two parts: a macroscopic (order N^) contribution and 
fluctuations around it (order iV°). Although the fluctu- 
ation component would be masked by the macroscopic 



7 



contribution, this latter is only non-zero in the superra- 
diant phase, such that the fluctuation-component could 
be experimentally accessed in the normal phase. 

We have identified three critical values Ai,A2,A3 for 
the coupling constant A. The value A2, Eq. (14), marks 



the phase transition, the area A > A2 has a non- vanishing 
macroscopic mode occupation |a|,|/?|, Eq. (13). Both 
the value of A 2 and of |q;| , |/3| coincide with results in 
the literatur(P2£l22l. The expression for Ai and A3 also 
agree with p revious calculations based an equations of 
motionglSESl^ however, only up to the first order in T 
(which is a small correction anyway). For example, our 
resuHP for Ai,A3, Eq. ([l8|(|l9]) differs from that of 
Dimeil^alet al. by a factor 2 in the term of ©(E^). The 
mean photon occupation number calculated for the nor- 
mal phase by Oztop^'^ et al. coincides with our result 
up to the first order in F. We ascribe these small dif- 
ference as to the different diagonalisation procedures for 
the effective Hamiltoni ans. 

As in previous work^Slifilsij^ ^j^g status of the region 
Al < A < A3, around the critical point making the 
normal-superradiant transition remains an open issue. In 
our approach, the effective Hamiltonian is not stable in 
that region. This could indicate a limitation of the simple 
mean field approach close to A2 in the dissipation model. 
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Appendix B: Master equation 



The three equations for the macroscopic occupation 



parameters a,a*,/3, Eq. (13), are 



Appendix A: Diagonalisation parameter 



In Eq. ( 15 1, the parameters for the diagonalisation are 

cos(7)(a; + e_) - cos(7)(w-£_) 



A = 
G = 

M = 

G2 = 



2Vw^/fr 

sin(7)(^ + g+). 



B 
D 



sin(7)(cj-£+)^ 



sin(7)(r2g + e_) ^ — _ sin(7)(-0g + £_) ^ 

cos(7)(f7o -e+) 



2^/rk>y/IZ 
cos(7)(rio + e+) 



2x/n^v^+ 



B2 
D2 



2Vr!oVeT 



(Al) 



In the normal phase, u) = VLq. The eigenenergies are 



2e 



2 -I- Wg ± \J (wq - uj'^ f + IGA^wwo, Wo > w 

2 + Wg =F \J (t^o ^ ^"^ f + ISA^WWg , Wg < W . 

(A2) 



= -2AW^y^+V^(w-z^) 



= -2A, / 1^ + Va^(^ + 



^ + «2 



= Aw|(^/S + ^/^) [l-^] -^^VP- 



Insertion of (B1),(B2) in (B3) produces 



= A 



2A./| 2A^A 
2j \ cj + ilj uj - if 



■ Wo 



VP- 



(Bl) 



(B2) 
(B3) 



k. 



(B4) 
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Using 1 — f = this expression can be simplified to 



= 



2(J - /?) 



V 



so that we get two solutions, /3 = or 



/3 = J 1 



r Wg WCJq 

16A2w " 4A2" 



J- 1 



/3, 
(B5) 

(B6) 



Using these values of a, we insert them into the pa- 

(2) 



rameters of Hg , Eq. (10), 
A = 



r2 + 40.2 



M = 



2V2v/w (16A2cj + (r2 + 4a;2) wq) 
8A2cj Wq 

r2 + 4cj2 + Y' 

Wo 96A4w2 - 2A2w (r2 + 4^2) 
'y ^ (r2 + 4a;2) (r2wo + 4w (4A2 + wwq)) 



(B7) 



Finally, the set of coupled differential equations for the 
P-representation, Eq. (27), reads (functions depends on 
t and x) 

a = ^W{cibi - di) - T{c2b2 - ^2) - 2{U + V) 

+ r(l-e*x)(B2 + £'^), 
bi=U -bi-W -bi-di, 

62 = t/ • 62 - T • 62 ■ ^2 , 
ci — U ■ ci — W • ci • di , 

C2 = V ■ C2 - f ■ C2 ■ d2 , 

di = 2U ■ di~W ■ dl 

+ r(l-e*x)(A2 + B2)^ 

d2^2V ■d2-T-dl 

+ r(l-e*x)(G2+L'2). (B8) 

The solution of this system in the long time limit is 

U + v/t/2 - (A2 + 52) (_i + g»x) p^r_ 



rfi(x) 



d2(x) 

^'i(x) 

a(i — > 00, x) 



V + v/V'^ - (-1 + e'^) + G2) TV 



T 



biix) = hix) = ci(x) = C2(x) = 0; 
2' 

f y/A'^ + - 2^252 (_i + 2e2»x) 

f Vi?^ - 21)2 (-1 + 2e2^x) 6-2 + G^^ R . 

(B9) 



First, we find a steady state solutions of the last six equa- 
tions, where time derivatives oi bi, Ci, di i G 1,2 are zero. 
Then we insert them into the differential equation for 
a{x,t), Eq. (B8), integrate, drop terms that are not 
dominant in case of t — >■ cx) and obtain the equation for 
a{t 00, x). 



Appendix C: System: Average in P-representation 

In order to obtain (^a\ai^, we have to use the rela- 
tion between the old and diagonal basis, because the P- 
representation has been evaluated in the diagonal basis. 



a\ai) = ( {Bdi + Ad\ + Dd2 + Gdl)x 

{Adi + Bd\ + Gd2 + Ddl)) . (CI) 



We have computed the P-representation in the interac- 
tion picture using the RWA, within which we obtain 



a\ai 



Tr J {{A^ + d\di + (G' + D^) d\d2 
l7i)(7i|® |72)(72| djfd-fi 



B^ + 



Jdjldji(^{A^ + B')jlj,-P 



(G^ + D^) 72*72 • P 



(C2) 



The P-representation is normalized. With our ansatz 



for P Eq. (28 1 we can represent the terms like 747* P as 



derivation with respect to di such that "fi"/iP = S^.P and 
do the integration J Pd'y^d'y'^ using symbolic numerical 
algebra to obtain. 



\^^^_D2 , n2 , /^2 , Tj2\bli't)ci{t)+di{t) 



\aA {t) = B'^ + + {A' + B^ 



di{tf 



(G^ + D^) _^2Wc2(0+c«2(t) 



rf2(t)2 



(C3) 



In a similar way. 



4^2) it) ^ Bl + Dl + (Al + P2) ^iW'^^W+'^iW 



/^2 ^ n2^ b2{t)c2{t) + d2{t) 



(C4) 
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H - 
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